clear
clear matrix
clear mata
set matsize 4000
set maxvar 10000

*Note: first change working directory replication package folder

run "scripts/programs/program_distlag.do"

/*********************************************************************************/
*Set program parameters
local radius = 20
local angle = 15
local length = 1
/**********************************************************************************/

/************************************** HU Data ***********************************/

use "data/`radius'-by-`length'km_`angle'deg/data_centroid.dta", clear
merge m:1 firenum using ///
	"data/`radius'-by-`length'km_`angle'deg/sample_assessors.dta", ///
	nogen

/************************************ Figure part 1 ***********************************/
/*********************************** Assessors Data ***********************************/


local controls lnt lnint majorroad tmissing wind_diff
local lags any_res count_res avgval_res val_res habpolytot wsimptot campground wilderness
local costlags tritot roadpcttot

label var any_res "No. res. props. > 0"
label var roadpcttot "Area < 0.5 km from road"
label var avgval_res "Avg. value res. props."
label var val_res "Total value res. props."
label var roadpcttot "Area < 0.5 km from road"
label var tmissing "$\Delta T$ missing"
label var lnt "Ln($\Delta T$ + 1)"
label var habpolytot "TES habitat"
label var wsimptot "Watershed importance"


foreach var in `lags' `costlags' {
	local j = `j' + 1
	local lab: variable label `var'
	local labels `"`labels' `j' "`lab'""'
	}

/************************************** Run regressions ***********************************/

local cond `"complexfire != 1 & year > 2011 & sample == 1"'

preserve
distlag y `controls'  if `cond', ///
	method(unrestricted) params(5) lagvars(`lags') costlags(`costlags') ///
	burnt(centroid) tmethod(log) fe(firenum) cluster(firenum)

gen spec = 1
gen data = 1
tempfile results
save `results', replace
restore

distlag y `controls'  if `cond', ///
	method(almon) params(2,5) lagvars(`lags') costlags(`costlags') ///
	burnt(centroid) tmethod(log) fe(firenum) cluster(firenum)

gen spec = 2
gen data = 1
append using `results'
save `results', replace

/************************************** MAKE FIGURE ***********************************/

// local legend = `"legend(order( `"Unrestricted weights"' 6 `"Almon weights, 5 km"' 7 `"Almon weights, 10 km"' 8 `"Linearly declining weights, 2 km"') size(vsmall))"'
local legend = `"legend(order(3 `"Unrestricted weights"' 4 `"Almon weights, quadratic"') size(small))"'
local axistitles `"ytitle(`"Change in probability of extinction"', margin(0 0 0 0) size(small)) xtitle(`"Lead distance (km)"', margin(0 0 2 1) size(small)) ylabel(,angle(0))"'
local graphlook `"graphregion(fcolor(white) margin(2 2 6 3)) `legend' yline(0, lpattern(dash) lcolor(red))"'

local color1 = "gs12"
local color2 = "navy"

local symbol1 = "Oh"
local symbol2 = "Th"

capture drop lag2
gen lag2 = lag + (spec - 1)*(.10)

capture drop panel
egen panel = group(data number)
capture label def varlabs `labels'
capture label values panel varlabs

forvalues i = 1/2 {

	if `i' == 1 local style `"`axistitles' `graphlook' "'

	local b`i' = `"(rspike ll ul lag2 if spec == `i' , lwidth(*1.5) color(`color`i''))"'
	local s`i' = `"(scatter beta lag2  if spec == `i', mlwidth(*1.5) mcolor(`color`i'') msymbol(`symbol`i'') lcolor(`color`i'') xlabel(,valuelabel) `style')"'
	local bars `bars' `b`i''
	local scatter `scatter' `s`i''
	}

twoway `bars' `scatter', by(panel, rescale yrescale cols(3) note("") colfirst iscale(0.5) graphregion(fcolor(white))) ///
	ysize(3) xsize(5) subtitle(,fcolor(gs14) size(medsmall)) scale(1.2)

graph export "Results/distlag_assessors.pdf", replace


/************************************ Figure part 2 ***********************************/
/************************************* Census Data ************************************/

use "data/`radius'-by-`length'km_`angle'deg/data_centroid.dta", clear
merge m:1 firenum using ///
	"data/`radius'-by-`length'km_`angle'deg/sample_census.dta", ///
	nogen

local controls lnt lnint majorroad tmissing wind_diff
local lags anyhu hu inc_pcincome inc_huxpcinc habpolytot wsimptot campground wilderness
local costlags tritot roadpcttot

label var anyhu "No. HU > 0"
label var hu "No. HU"
label var inc_pcincome "Per cap. income"
label var inc_huxpcinc "No. HU x per cap. inc."
label var roadpcttot "Area < 0.5 km from road"
label var tmissing "Delta T missing"
label var lnt "Ln($\Delta T$ + 1)"
label var habpolytot "TES habitat"
label var wsimptot "Watershed importance"

local labels
local j
foreach var in `lags' `costlags' {
	local j = `j' + 1
	local lab: variable label `var'
	local labels `"`labels' `j' "`lab'""'
	}

/************************************** Run regressions ***********************************/

local cond `"complexfire != 1 & sample == 1"'

preserve
distlag y `controls'  if `cond', ///
	method(unrestricted) params(5) lagvars(`lags') costlags(`costlags') ///
	burnt(centroid) tmethod(log) fe(firenum) cluster(firenum)

gen spec = 1
gen data = 1
tempfile results
save `results', replace
restore

distlag y `controls'  if `cond', ///
	method(almon) params(2,5) lagvars(`lags') costlags(`costlags') ///
	burnt(centroid) tmethod(log) fe(firenum) cluster(firenum)

gen spec = 2
gen data = 1
append using `results'
save `results', replace

/************************************** MAKE FIGURE ***********************************/

/// local legend = `"legend(order( `"Unrestricted weights"' 6 `"Almon weights, 5 km"' 7 `"Almon weights, 10 km"' 8 `"Linearly declining weights, 2 km"') size(vsmall))"'
local legend = `"legend(order(3 `"Unrestricted weights"' 4 `"Almon weights, quadratic"') size(small))"'
local axistitles `"ytitle(`"Change in probability of extinction"', margin(0 0 0 0) size(small)) xtitle(`"Lead distance (km)"', margin(0 0 2 1) size(small)) ylabel(,angle(0))"'
local graphlook `"graphregion(fcolor(white) margin(2 2 6 3)) `legend' yline(0, lpattern(dash) lcolor(red))"'

local color1 = "gs12"
local color2 = "navy"
local color3 = "dkgreen"
local color4 = "dkorange"

local symbol1 = "Oh"
local symbol2 = "Th"
local symbol3 = "Sh"
local symbol4 = "Dh"

capture drop lag2
gen lag2 = lag + (spec - 1)*(.10)

capture drop panel
egen panel = group(data number)
capture label def varlabs `labels'
capture label values panel varlabs

forvalues i = 1/2 {

	if `i' == 1 local style `"`axistitles' `graphlook' "'

	local b`i' = `"(rspike ll ul lag2 if spec == `i' , lwidth(*1.5) color(`color`i''))"'
	local s`i' = `"(scatter beta lag2  if spec == `i', mlwidth(*1.5) mcolor(`color`i'') msymbol(`symbol`i'') lcolor(`color`i'') xlabel(,valuelabel) `style')"'
	local bars `bars' `b`i''
	local scatter `scatter' `s`i''
	}

twoway `bars' `scatter', by(panel, rescale yrescale cols(3) note("") colfirst iscale(0.5) graphregion(fcolor(white))) ///
	ysize(3) xsize(5) subtitle(,fcolor(gs14) size(medsmall)) scale(1.2)

graph export "Results/distlag_census.pdf", replace
